1/1 


* 


D-A184  454 


UNCLASSIFIED 


SINILARITV  SOLUTION  FOR  PERTURBATION  OF  A  SRRONG  BLAST 
NAVE  WITH  RAD  I AT  1 0  (U>  MISSION  RESEARCH  CORP  SANTA 

BARBARA  CA  0  HCCARTOR  ET  AL  21  DEC  BE  HRC-R-767 
DNA-TR-87-28  DNA881 -85-C-8825  F/C  19/11 


AD-A184  454 


DIB  fiii  Cutur  ' 

DNA-TR-87-38 


SIMILARITY  SOLUTION  FOR  PERTURBATION  OF  A 
STRONG  BLAST  WAVE  WITH  RADIATION  DIFFUSION  IN 
A  STRATIFIED  MEDIUM 


G.  McCartor 
W.  R.  Wortman 

Mission  Research  Corporation 

P.  0.  Drawer  719 

Santa  Barbara,  CA  93102-0719 


31  December  1986 


Technical  Report 


DTIC 


ELECTE 
SEP  1  0  1987 


CONTRACT  No.  DNA  001-85-C-0035 


Approved  for  public  release; 
distribution  is  unlimited. 


THIS  WORK  WAS  SPONSORED  BY  THE  DEFENSE  NUCLEAR  AGENCY 
UNDER  RDT&E  RMC  CODE  B3220857602  RZ  RA  00001  25904D. 


Prepared  for 
Director 

DEFENSE  NUCLEAR  AGENCY 
Washington,  DC  20305-1000 


Destroy  this  report  when  it  is  no  longer  needed.  Do  not  return 
to  sender. 

PLEASE  NOTIFY  THE  DEFENSE  NUCLEAR  AGENCY 
ATTN:  TITL,  WASHINGTON,  DC  20305  1000,  IF  YOUR 
ADDRESS  IS  INCORRECT,  IF  YOU  WISH  IT  DELETED 
FROM  THE  DISTRIBUTION  LIST,  OR  IF  THE  ADDRESSEE 
IS  NO  LONGER  EMPLOYED  BY  YOUR  ORGANIZATION. 


CUT  HERE  AND  RETURN 


DISTRIBUTION  LIST  UPDATE 


This  mailer  is  provided  to  enable  DNA  to  maintain  current  distribution  lists  for  reports.  We  would 
appreciate  your  providing  the  requested  information. 

□  Add  the  individual  listed  to  your  distribution  list. 

□  Delete  the  cited  organization/individual. 

□  Change  of  address. 

NAME:  _ 

ORGANIZATION:  _ 

OLD  ADDRESS  CURRENT  ADDRESS 


TELEPHONE  NUMBER:  S _ L_ 

SUBJECT  AREA(s)  OF  INTEREST: 


DNA  OR  OTHER  GOVERNMENT  CONTRACT  NUMBER: _ 

CERTIFICATION  OF  NEED-TO-KNOW  BY  GOVERNMENT  SPONSOR  (if  other  than  DNA): 

SPONSORING  ORGANIZATION:  _ 

CONTRACTING  OFFICER  OR  REPRESENTATIVE: _ 

SIGNATURE: _ 


Director 

Defense  Nuclear  Agency 
ATTN:  mm  TITL 
Washington,  DC  20305-1000 


Director 

Defense  Nuclear  Agency 
ATTN:  H  TITL 
Washington,  DC  20305-1000 


UNCLASSIFIED 


REPORT  DOCUMENTATION  PAGE 


Form  Approved 
OMB  No.  0704-0188 
Exp.  Date:  Jun 30,  1988 


la.  REPORT  SECURITY  CLASSIFICATION 

UNCLASSIFIED 

2a.  SECURITY  CLASSIFICATION  AUTHORITY 

N/A  since  Unclassified 

2b.  DECLASSIFICATION /DOWNGRADING  SCHEDULE  ll 

N/A  since  Unclassified 

_ 1 

4  PERFORMING  ORGANIZATION  REPORT  NUMBER(S)  | 

MRC-R-767 

6a.  NAME  OF  PERFORMING  ORGANIZATION 

6b.  OFFICE  SYMBOL 

Mission  Research  Corporation 

(If  applicable) 

6c.  ADDRESS  (Gfy,  State,  and  ZIP  Code) 

P.  0.  Drawer  719 

Santa  Barbara,  CA  93102-0719 

8a.  NAME  OF  FUNDING  /  SPONSORING 

8b.  OFFICE  SYMBOL 

ORGANIZATION 

(If  applicable) 

RAAE/Schwartz 

8c.  ADDRESS  (Gfy,  State,  and  ZIP  Code) 

1b.  RESTRICTIVE  MARKINGS 


3  DISTRIBUTION /AVAILABILITY  OF  REPORT 

Approved  for  public  release; 
distribution  is  unlimited. 


5.  MONITORING  ORGANIZATION  REPORT  NUMBER(S) 

DNA-TR-87-38 


7a.  NAME  OF  MONITORING  ORGANIZATION 

Director 

Defense  Nuclear  Aqenc 


7b  ADDRESS  (City,  State,  and  ZIP  Code) 

Washington,  DC  20305-1000 


9.  PROCUREMENT  INSTRUMENT  IDENTIFICATION  NUMBER 

DNA  001-85-C-0035 


10  SOURCE  OF  FUNDING  NUMBERS 


PROGRAM  PROJECT 

ELEMENT  NO  NO 

62715H  RZ 


1 1  TITLE  (Include  Security  Classification) 

SIMILARITY  SOLUTION  FOR  PERTURBATION  OF  A  STRONG  BLAST  WAVE  WITH  RADIATION 
DIFFUSION  IN  A  STRATIFIED  MEDIUM 


12.  PERSONAL  AUTHOR(S) 

McCartor,  G.;  Wortman,  W.  R. 


13a  TYPE  OF  REPORT  13b  TIME  COVERED  14  DATE  OF  REPORT  (Year,  Month,  Day)  115  PAGE  COUNT 

Technical  from  860712  TQ  861231  861231  50 


16  SUPPLEMENTARY  NOTATION 

This  work  was  sponsored  by  the  Defense  Nuclear  Agency  under  RDT&E  RMC  Code 
B3220857602  RZ  RA  00001  25904D. 


17  COSATI  CODES  18  SUBJECT  TERMS  (Continue  on  reverse  if  necessary  and  identify  by  block  number) 

field  group  sub-group  B1  astwaves  Flow  Fields 

19  9  Radiation  Hydrodynamics 

__  _  (_ow  AHitucje  Bursts 


ABSTRACT  ( Continue  on  reverse  if  necessary  and  identify  by  block  number ) 


i 


similarity  solution  for  the  first  order  perturbation  of  a  nearly  spherical  strong  blast 
wave  generated  by  a  point  source  in  a  mildly  stratified  medium  is  found.  The  hydrodynamic 
equations  include  a  radiation  diffusion  approximation  to  radiation  transport  which  is  of  such 
a  character  to  allow  both  the  unperturbed  and  perturbed  solutions  to  exhibit  a  similarity 
form.  The  cyl indrical ly  symmetric  solution  provides  the  flow  in  both  radial  and  polar 
directions.  Implications  for  circulation  in  flow  fields  of  low-altitude  atmospheric 
explosions  are  drawn. 


20  DISTRIBUTION /AVAILABILITY  OF  ABSTRACT 
□  UNCLASSIFIFD/UNLIMITED  03  SAME  AS  RPT  □  DTiC  USERS 


22a  NAME  OF  RESPONSIBLE  INDIVIDUAL 

Sandra  E.  Youn 


OO  FORM  1473, 84  MAR  83  APR  edition  may  be  used  until  exhausted 

All  other  editions  are  obsolete. 


21  ABSTRACT  SECURITY  CLASSIFICATION 

UNCLASSIFIED 


22b  TELEPHONE  (Include  Area  Code)  22c  OFFICE  SYMBOL 

(202)  325-7042  DNA/CSTI 


il  exhausted  SECURITY  CLASSIFICATION  OF  THIS  PAGE 

UNCLASSIFIED 


TABLE  OF  CONTENTS 


Section  Page 

LIST  OF  ILLUSTRATIONS  iv 

1  INTRODUCTION  1 

2  RADIATION  HYDRODYNAMIC  EQUATIONS  4 

2.1  General  Case  4 

2.2  Zeroth  Order  5 

2.3  First  Order  10 

3  SOLUTIONS  TO  SIMILARITY  EQUATIONS  18 

4  VORTICITY  AND  VELOCITY  AT  PRESSURE  EQUILIBRIUM  23 

5  DISCUSSION  28 

6  LIST  OF  REFERENCES  30 

Appendix 

EXTENSION  OF  SIMILARITY  PERTURBATION  TECHNIOUE  31 

A. 1  GENERAL  FORM  31 

A. 2  SPECIFIC  EQUATIONS  AT  2ND  ORDER  37 


LIST  OF  ILLUSTRATIONS 


Figure  Page 

1  Density  similarity  profile.  8 

2  Velocity  similarity  profile.  8 

3  Pressure  similarity  profile.  9 

4  Temperature  similarity  profile.  9 

5  Perturbed  radial  velocity  similarity  profile.  21 

6  Perturbed  polar  velocity  similarity  profile.  21 

7  Perturbed  temperature  similarity  profile.  22 

8  Perturbed  density  similarity  profile.  22 

9  Vertical  velocity  versus  radius  at  pressure  equilibrium  26 

time  at  e=0  (vertical  direction). 

10  Vertical  velocity  versus  radius  at  pressure  equilibrium  26 

time  at  e=w/2  (horizontal  direction). 

11  Vertical  velocity  versus  radius  at  pressure  equilibrium  27 

time  at  0=ir/4  (45*  off  vertical  direction). 


SECTION  1 
INTRODUCTION 


The  distortion  of  a  strong  blast  wave  propagating  in  an  inhomo¬ 
geneous  medium  has  been  the  subject  of  a  variety  of  investigations.  A 
summary  of  the  older  work  is  given  by  Zel'dovich  and  Raizer.1  Kompaneets2 
looked  at  the  problem  of  the  propagation  of  an  initially  spherical  blast 
wave  from  a  point  source  in  an  atmosphere  stratified  in  one  dimension  with 
exponentially  changing  density.  He  assumed  strictly  radial  flow  along 
with  a  hypothesis  of  uniform  pressure  in  the  interior  of  the  expanding 
region  driving  a  thin  shock  front  containing  all  the  swept  up  mass.  This 
resulted  in  a  simple  analytic  solution  for  the  shock  front  motion. 
Raizer3  found  a  self-similar  solution  for  an  impulsive  source  for  planar 
shock  propagation  in  an  exponential  medium;  Grover  and  Hardy4  have 
reported  this  solution  along  with  numerical  results  including  spherical 
geometry  (with  the  source  at  the  center  of  the  spherical  density  varia¬ 
tion).  Laumbach  and  Probstein5  have  given  an  improved  version  of 
Kompaneet's  calculation  by  avoiding  any  ad  hoc  assumption  about  the  pres¬ 
sure  field  while  retaining  the  assumption  of  radial  flow.  They  do  assume 
a  strong  blast  wave  and  solve  the  hydrodynamic  equations  for  a  point 
source  under  an  assumption  which  takes  the  bulk  of  the  mass,  from  inside 
the  current  polar  angle,  to  be  concentrated  near  the  shock.  This  is  used 
to  truncate  expansions  behind  the  shock.  The  result  is  then  an  approxima¬ 
tion  to  both  the  shock  location  and  profiles  of  hydrodynamic  variables 
behind  the  shock  which  is  a  function  of  polar  angle,  scale  height  and 
energy. 


In  the  limit  of  a  homogeneous  atmosphere  (horizontal  propaga¬ 
tion),  the  solution  of  Laumbach  and  Probstein  provides  a  good  approxima¬ 
tion  to  the  well-known  strong  blast  wave  solution  given  by  Taylor6  and 
Sedov.7  The  Taylor-Sedov  solution  can  be  put  in  similarity  form  so  that 
the  hydrodynamic  variables  at  any  time  and  location  can  be  expressed  as  a 
product  of  a  power  to  time  and  a  fixed  profile  which  is  a  function  of  only 
a  similarity  variable  which  may  be  taken  as  the  radius  over  the  radius  to 
the  shock.  Obviously  the  resulting  separation  of  variables  greatly  simp¬ 
lifies  the  solution  of  the  problem.  However,  as  indicated  by  Laumbach  and 
Probstein,  the  general  case  of  an  exponentially  changing  ambient  density 
does  not  admit  to  a  similarity  solution. 

The  Taylor-Sedov  solution  is  not  fully  useful  for  the  realistic 
problem  of  a  point  energy  source  in  a  uniform  atmosphere  because  it 
assumes  that  radiation  energy  and  transport  are  insignificant.  Elliott8 
has  studied  the  circumstances  under  which  this  restriction  may  be  lifted. 
He  finds  that  when  radiation  transport  is  added  to  the  energy  conservation 
equation  by  the  radiation  diffusion  approximation,  the  resulting  equations 
sometimes  still  allow  a  similarity  solution.  If  the  temperatures  are  high 
enough  for  the  radiant  flux  to  affect  significantly  the  hydrodynamics  but 
sufficiently  low  that  the  radiation  energy  and  pressure  are  small  compared 
to  the  material  values,  a  particular  form  for  the  diffusion  coefficient 
will  give  a  similarity  solution.  Remarkably  enough,  the  form  needed  is 
quite  similar  to  that  seen  experimentally.  We  shall  exploit  this  fact  in 
constructing  approximate  solutions  for  flow  fields  for  a  class  of  atmo¬ 
spheric  explosions. 

In  this  paper  the  flow  fields  influenced  by  radiation  diffusion 
behind  a  strong  shock  from  an  energy  point  source  will  be  found  for  an 
atmosphere  mildly  stratified  in  density.  That  is,  the  first  order  pertur¬ 
bation  in  the  flow  induced  by  a  noninfinite  scale  height  will  be  deter¬ 
mined.  It  is  found  that  both  the  uniform  atmosphere  case,  done  by 


Elliott8,  and  the  first  order  perturbation,  with  an  expansion  coefficient 
inversely  proportional  to  the  scale  height,  allow  similarity  solutions  for 
which  the  dependent  variables  can  be  expressed  as  the  product  of  a  power 
of  time,  functions  of  the  similarity  variable  and  simple  angular  func¬ 
tions.  In  fact,  it  can  be  shown  that  such  similarity  solutions  exist  to 
all  orders  of  the  expansion  although  this  paper  will  be  concerned  only 
with  the  first  order.  However,  each  order  requires  a  different  similarity 
transformation. 


The  motivation  for  pursuit  of  this  problem  lies  in  the  determin¬ 
ation  of  the  residual  flow  fields  which  influence  the  rise  and  distortion 
of  fireballs  from  low-altitude  explosions.  A  realistic  calculation  of  all 
features  of  such  flow  requires  much  more  than  the  current  calculation.  At 
very  early  times,  explosion  dynamics  are  dominated  by  radiation  transport 
which  is  not  diffusive  and  for  which  the  radiant  energy  is  important. 
Once  the  radiation  sphere  cools  sufficently,  a  shock  is  launched  from  its 
surface  of  the  radiation  sphere.  This  strong  shock  propagates  through  the 
ambient  atmosphere  in  a  manner  described  by  the  current  calculation.  In  a 
stratified  atmosphere  the  propagation  is  distorted  by  a  more  rapid  motion 
into  lower  densities.  Ultimately,  the  shock  propagation  is  weakened  by 
spherical  divergence  so  that  it  is  no  longer  a  strong  shock.  At  this 
point  the  current  calculation  will  no  longer  be  appropriate.  Ultimately 
the  expansion  of  this  region  will  cease  when  the  pressure  is  reduced  to 
ambient  levels.  When  this  pressure  equilibrium  fireball  condition  is 
reached,  there  will  be  a  portion  of  the  fluid  through  which  the  strong 
shock  propaqated  which  later  is  nearly  isothermal  due  to  diffusion  of 
radiation  behind  the  shock.  It  is  precisely  this  region  for  which  the 
circulation,  induced  by  mild  atmospheric  stratification,  which  can  be  cal¬ 
culated  using  the  methods  of  this  paper.  This  region  will  have  its  circu¬ 
lation  frozen  and  convected  with  the  flow.  This  circulation  is  important 
in  providing  the  character istic  torus  at  later  times  as  the  now  underdense 
pressure  equilibrium  fireball  undergoes  buoyant  rise.  For  the  current 
calculation,  only  the  strong  shock  flow  fields  will  be  found. 


SECTION  2 


RADIATION  HYDRODYNAMIC  EQUATIONS 


GENERAL  CASE. 


The  equations  for  inviscid  hydrodynamic  flow  for  an  ideal  fluid 
in  the  presence  of  a  radiant  energy  flux  are  taken  to  be 


3p/3t  +  v.(pv)  =  0 
3 v/st  +  (v*v)v  =  -  (vp)/p 
3p/3t  +  ypV  •  V  +  V»Vp  +  (y-l)tf*F 


where  y  is  the  adiabatic  exponent  and  F  is  the  radiant  flux.  In  the  radi¬ 
ation  diffusion  approximation,  the  flux  is  proportional  to  the  temperature 
gradient  as 


F  =  -  n  v  T 


where  n  is  the  diffusion  coefficient  which  is  inversely  proportional  to 
an  opacity,  K,  given  by 


K  =  Ifiall 

3pn 


where  a  is  the  Stefan-Boltzmann  constant.  This  radiation  diffusion 


approximation  is  valid  when  the  photon  mean  free  paths  are  small  compared 
with  the  hydrodynamic  dimensions  of  interest  in  the  problem  and  the 
radiant  energy  and  pressure  are  small  compared  with  the  corresponding 
hydrodynamic  quantities.  The  latter  conditions  permit  the  introduction  of 
an  ideal  fluid  equation  of  state  which  completes  the  set  of  equations. 
Apart  from  early  times,  these  equations  should  provide  a  reasonable 
description  for  low-altitude  explosions. 


The  above  equations  appear  everywhere  except  at  the  shock  front. 
At  the  shock  the  hydrodynamic  variables  must  meet  Hugoniot  jump  condi¬ 
tions9  which  enforce  mass,  momentum  and  energy  conservation  through  the 
shock.  Under  the  conditions  assumed,  these  will  now  be  expressed  in 
terms  of  the  local  ambient  density  p,  and  the  local  shock  velocity  vc  as 

d  S 


P  (  V 
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V 
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=  0 


p  =  pa  V-*. 


yP 

(y-l)p 


+ 


(6) 

(7) 

(8) 

(9) 


where  it  is  assumed  that  F  =p=v  =v=0  outside  the  shock.  Here  p=0 
is  the  strong  shock  assumption  while  F=0  indicates  that  no  radiant  energy 
escapes  the  region  so  that  total  energy  inside  the  shock  will  be  conserv¬ 
ed.  Note  that  ii  and  4-  indicate  parallel  and  perpendicular  to  the  local 
shock  surface  so  that  the  location  of  the  shock  front,  including  it., 
velocity,  is  a  part  of  the  solution  to  the  problem.  In  principle,  the 
above  equations  will  provide  for  strong  shock  motion  generated  for  any 
ambient  density  p  and  opacity  once  the  condition  of  total  energy  conser- 
vation  is  imposed.  However,  solutions  are  available  only  for  a  limited 
set  of  cases  where  similarity  conditions  are  met.  (For  the  case  where  the 
strong  shock  condition  is  dropped,  numerical  solutions  have  been  found  by 
Brode10  for  the  case  of  uniform  ambient  density  and  no  diffusion.) 


2.2  ZEROTH  ORDER. 


In  the  limit  of  uniform  ambient  density,  which  is  the  zeroth 
order  case  for  an  expansion  in  a  parameter  inversely  proportional  to  the 
scale  height,  the  solution  for  spherically  symmetric  strong  shock  problem 


is  known  from  Taylor6  and  Sedov7  for  no  diffusion.  It  is  evident  from 
dimensional  arguments  that  there  is  a  similarity  solution  which  has  the 
properties  that  the  shock  radius  R(t)  is  proportional  to  t2/5  ,  while  the 
density,  fluid  velocity,  and  pressure  can  be  expressed  as  a  power  of  t 
(the  powers  are  0,  -3/5,  and  -6/5,  respectively)  multiplied  by  profiles 
which  are  functions  only  of  a  similarity  variable  which  can  be  taken  as 
the  ratio  of  the  radius,  r)  to  the  shock  radius,  R(t).  Relatively  simple 
closed  form  solutions  for  the  profiles  exist. 

When  radiation  diffusion  is  included  for  the  constant  ambient 
density  case,  it  is  not  generally  possible  to  find  similarity  solutions 
for  the  equations.  However,  Elliott8  has  demonstrated  that  for  a  diffu¬ 
sion  coefficient  with  a  suitable  dependence  on  temperature,  the  similarity 
form  can  be  maintained.  It  is  easy  to  show  that  the  diffusion  term  in  the 
energy  equation  (3)  will  give  the  same  separation  into  a  power  of  time 
times  a  scaled  profile,  as  do  the  other  terms,  using  the  same  powers  as 
indicated  by  Taylor  provided  that  the  diffusion  coefficient,  n,  is  propor¬ 
tional  to  1  /6th  power  of  the  temperature,  T.  It  may  be  any  function  of 
the  density.  This  means  the  opacity  must  be  proportional  to  T17/6  but  may 
have  any  density  dependence.  It  turns  out  that  for  temperatures  below 
about  10  eV,  the  range  of  interest  for  the  current  problem,  the  Rosseland 
mean  opacity  rises  sharply  approximately  like  T3  while  decreasing  with 
density.  For  calcul ational  purposes,  the  diffusion  coefficient  will  be 
taken  as 

n  =  a  pk  T1/6  (10) 

where  k=  -1.6  as  estimated  for  air  densities  between  10"3  and  10-4 
gm/cm3  which  gives  a  «  8.7  x  10'5  in  cgs  units. 

Elliott  has  given  a  technique  for  solving  the  similarity  equa¬ 
tions  which  involves  carrying  out  the  first  integral  of  the  energy  equa¬ 
tion  analytically  and  making  a  change  of  variables  which  converts  all  the 


unknown  boundary  values  to  appear  at  the  (unknown)  shock  location.  The 
equations  are  then  integrated  outward  from  the  origin  to  find  the  scaled 
location  of  the  shock  front  at  which  the  Hugoniot  conditions  are  met.  The 
only  drawback  to  this  procedure  is  that  the  required  scaling  of  the 
diffusion  coefficient  is  dependent  on  the  solution  so  it  is  necessary  to 
first  find  a  solution  with  a  scaled  coefficient  and  then  determine  the 
physical  coefficient  to  which  it  corresponds.  As  a  result  it  necessary  to 
iterate  to  find  the  solution  for  a  given  physical  coefficient.  In  prac¬ 
tice,  this  is  not  a  problem  since  there  is  a  smooth  relation  between  the 
scaled  and  physical  coefficients.  A  part  of  the  solution  also  consists  of 
finding  the  constant  of  proportionality  between  the  shock  radius,  R,  and 
t2/5.  This  is  a  function  of  E  (total  energy),  p  ,  y,  and  n  and  it  can  be 

Q 

determined  by  integrating  the  solution  profiles  for  internal  and  kinetic 
energy  and  equating  the  sum  to  the  total  energy.  For  Taylor's  case  (n  = 
0)  it  is  known  that  for  y  =  1.2  (which  is  a  reasonable  value  for  ionized 
air)  the  energy  relation  is 

E  =  11.0  PaR2R3  (11) 
and 

R  =  a  t2/5  which  gives  a  in  terms  of  E. 

For  a  pa  of  5  •  10-l+  gm/cm3  (which  corresponds  to  approximately 
one  scale  height  above  sea  level)  and  an  energy  of  40  kT  (or  1.67  •  1021 
ergs)  the  (dimensionless)  coefficient  of  11.0  is  reduced  to  10.2  by  inclu¬ 
sion  of  diffusion  indicating  a  moderate  effect  near  the  shock  which  serves 
to  increase  the  shock  speed.  (The  example  used  is  similar  to  that  labeled 
as  K  =  0.5  in  Elliott's  case.)  Although  the  effects  of  diffusion  near  the 
shock  are  modest,  the  effects  in  the  center  of  the  fireball  are  pronounced 
as  the  temperature  is  forced  to  be  fairly  uniform  over  the  bulk  of  the 
volume.  This  is  illustrated  by  the  plots  of  density,  velocity,  pressure 
and  temperature  shown  in  Figures  1-4.  The  plots  show  the  profiles  taken 
as  a  function  of  the  similarity  variable  x  =  r/R(t)  where 


V, 

4 


p(r,t)  =  Pl (x)  (12) 
v(r,t)  =  v^x)  t'3/5  (13) 
p(r,t)  =  Pl(x)  f6/5  (14) 
T(r,t)  =  Tx(x)  f6/5  (15) 

where  T  is  not  independent  but  given  by  the  ideal  gas  equation  of  state 
p  =  pR*T  (16) 

★ 

Here  R  is  the  gas  constant  given  by  the  product  of  Boltzmann's  constant 
and  Avogadro's  number  over  the  mean  molecular  weight  of  the  constitutents 
which  has  been  taken  as  thirty  for  the  numerical  example. 

The  solutions  including  diffusion  differ  from  the  Taylor's  case 
in  that  the  temperature  now  is  nearly  uniform  in  the  interior  while  the 
pressure  and  density  go  to  nonzero  values  at  the  origin.  Near  the  origin 
the  pressure  goes  like  the  fourth  power  of  x  plus  a  constant,  the  density 
goes  like  the  second  power  of  x  plus  a  constant  while  the  velocity,  to 
lowest  order,  goes  like  the  third  power  of  x.  Near  the  shock,  the 
solutions  are  similar  to  those  without  diffusion  although  the  speed  of  the 
front  in  increased  somewhat  by  the  outward  flow  of  radiant  energy. 

From  the  present  considerations,  it  is  important  that  there 
exists  an  easily  obtained  solution  which  strongly  resembles  that  for  a 
physical  fireball  which  can  serve  a  the  basis  for  a  perturbation  expansion 
to  determine  the  character  of  the  effects  of  mild  stratification  of  the 
ambient  medium. 

2.3  FIRST  ORDER. 

In  ambient  density  which  is  not  constant,  the  solution  given 
above  is  no  longer  correct.  However,  if  the  density  changes  only  slightly 
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over  the  region  defined  by  the  shock  position,  it  is  expected  that  this 
solution  will  be  approximately  correct.  This  fact  can  be  exploited  for 
the  case  of  the  atmosphere,  whose  density  variation  is  almostly  exclusive¬ 
ly  in  the  vertical  direction  and  can  be  described  locally  as  exponential 
with  a  constant  characteristic  scale  height,  H,  which  is  imposed  by  an 
equilibrium  between  gravity  and  ambient  pressure  forces.  The  dynamical 
Equations  1-3  are  still  appropriate  for  this  nonconstant  density  provided 
that  gravitational  forces-  are  not  important  for  development  of  the  strong 
shock  and  associated  flow.  Since  strong  shock  conditions  are  imposed, 
this  is  true.  However,  the  effects  of  an  exponential  ambient  density  will 
alter  the  Hugoniot  conditions  in  a  mild  way  (so  long  as  R(t)  «  H)  exclu¬ 


sively  through  the  altitude  dependence  of  p  .  Since  the  density  is  not 

a 


spherically  symmetrical,  the  resulting  flow  will  now  have  a  polar  angle 
dependence,  including  a  nonradial  velocity,  and  the  possibility  of  vortic- 
ity  where  there  was  none  before. 


Previous  work1-5  indicates  that  when  a  strong  shock  propagates 
in  an  exponential  atmosphere,  the  motion  into  decreasing  density  is 
enhanced  while  that  in  increasing  density  is  diminished.  The  initially 
spherical  pulse  is  thus  distorted  by  the  upward  expansion  which  can,  for  a 
small  H,  result  in  the  upwardly  moving  shock  accelerating  all  the  way 
through  the  top  of  the  atmosphere  (where  the  hydrodynamic  equations  are, 
of  course,  not  valid).  For  a  large,  but  finite,  scale  height,  a  mild  dis¬ 
tortion  of  the  spherical  pulse  is  expected.  This  distortion  will  now  be 
described  as  a  small  perturbation  about  the  constant  density  solution  and 
an  expansion  for  perturbed  solution  will  be  found. 


The  ambient  density  and  its  leading  term  is 


p  =  p  exp  (-z/H)  =  Dfl  exp(-r  cose/H)  *  p  (1-r  cose/H)  (17) 


where  p  is  the  ambient  density  at  z=0,  9  is  the  polar  angle  and  r  is  the 


spherical  radius.  The  form  of  the  expanded  density  suggests  the  expansion 


of  the  perturbed  solution  in  a  Fourier  series  in  order  to  define  the  angu¬ 
lar  dependence.  When  this  is  done,  it  is  readily  found  that  the  series 
for  all  the  dependent  variables,  including  the  boundary  conditions,  trun¬ 
cate  to  the  following  form  provided  that  higher  order  terms  are  ignored: 


p  =  pi  +  f  (r,t)  cose  (18) 

p 

P  =  Pi  +  fp(r,t)  cose  (19) 

vr  =  vj  +  fy(r,t)  cose  (20) 

v9  ■  Vr’t)  sin9  (21) 


The  position  of  the  perturbed  shock  front  is  described  by  an  auxilliary 
function  G(r,e,t)=0  where 

G  =  R(t)  -  r  +  fR(t)  cose  (22) 


so  fR  is  the  perturbation  in  shock  front  location  in  the  vertical  direc¬ 
tion.  The  (linearized)  equations  which  the  f(r,t)  functions  obey  are: 
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Here  the  subscript,  one,  has  been  dropped  from  the  leading  order  terms. 
The  boundary  conditions  to  this  order,  as  applied  at  the  perturbed  shock 
location,  are 

v  f„  ■  Rf„  (27) 
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where  the  temperature  perturbation  fy  has  been  used  as  an  alternative  to 
the  pressure  through  the  first  order  terms  in  the  equation  of  state.  It 
must  be  pointed  out  that  higher  order  terms,  such  as  f2,  will  contain 
Fourier  contributions  which  are  independent  of  anqle  and  which  will  thus 
couple  the  zeroth  and  second  order  terms.  However,  when  only  first  order 
terms  are  kept,  the  only  terms  constant  in  angle  are  those  from  the  zeroth 
order. 


The  perturbation  equations  for  four  f  functions  consist  of  a 
linear  set  of  four  coupled  homogeneous  second  order  partial  differential 
equations  in  two  dependent  variables  subject  to  four  linear  but  inhomoge¬ 
neous  boundary  conditions  applied  at  a  position  fixed  by  fR  which  is  to  be 
determined.  The  task  of  solvinq  would  be  formidable  except  for  one  fact. 
It  is  possible  to  find  a  similarity  solution  to  the  perturbed  equations  in 
a  manner  completely  analogous  to  that  seen  at  zeroth  order!  It  is  a 
straightforward  matter  to  verify  that  when  the  f's  are  written  as 


f  (r,t)  =  t2/5  g  (x) 


fpfr.t)  =  t"4/5  gp(x) 


f  (r’t)  =  t"1/5  g  (x) 


f0(r,t)  =  t-1/5  q0(x) 


fR(t)  =  t4/5  gF 


where  x  =  r/R(t)  as  before,  the  equations  are  reduced  to  similarity  form 
in  the  variable  x.  Consequently,  the  problem  is  reduced  to  solving  a  set 
of  coupled  second  order  linear  differential  equations  in  one  dependent 
variable.  The  equations  are 
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where  the  equation  of  state  at  this  order  is 

P-  =  ^  .  (40) 
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The  Hugonoit  conditions  applied  at  x=i  are 
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In  the  presence  of  the  perturbation,  the  total  energy  must  be 
unchanged.  This  condition  is  automatically  imposed  at  this  order  by  the 
symmetry  of  the  angular  dependence.  The  total  momentum  must  also  be 
constant  (and  equal  to  zero)  in  the  presence  of  the  perturbation.  This 
condition  is  not  automatically  satisfied  and  it  serves  as  an  additional 
condition.  The  vertical  momentum  integral  is 

1 

PZ  “  9rv(1 )  p(1)  +  a  /  x2dx  [gpv  +  Pgy  -  2Pg0]  (45) 

The  second  order  energy  equation  can  be  trivially  changed  into 
two  first  order  equations  if  the  dependent  variable 

fq  =  df  T/ dx  (46) 

is  introduced  and  the  second  derivative  of  fy  is  replaced  by  dfR/dx. 

There  are  then  five  first  order  equations  with  two  unknown  parameters 
which  fix  the  values  at  x=l.  These  are  fR(l)  and  gR.  Note  that  qR  along 


with  the  zeroth  order  solutions  fully  fix  the  Hugoniot  conditions  at  this 
order.  In  principle,  there  can  be  a  family  of  solutions  in  these  two 
parameters  for  which  only  parameter  can  be  specified  independently  and 
still  allow  the  momentum  to  vanish.  Consequently,  one  more  condition  can 
be  imposed  on  physical  grounds.  The  required  condition  is  that  the 

perturbation  in  the  temperature,  fy,  vanish  at  the  origin.  If  it  does 
not,  there  will  be  a  discontinuity  in  temperature  there  in  crossing  the 
horizontal  plane  of  the  source  since  the  coefficient  of  fy  is  cose  which 
changes  sign  across  the  plane.  Now  it  is  possible,  in  principle,  to  solve 
the  equations  by  exploring  the  solutions  in  the  parameters  fq(l)  and  9R 
and  selecting  that  which  causes  the  total  momentum  and  fy(0)  to  vanish. 


SECTION  3 

SOLUTIONS  TO  SIMILARITY  EQUATIONS 


The  set  of  five  coupled  linear  first  order  differential  equa¬ 
tions,  which  contain  the  known  zeroth  order  functions  as  variable  coeffi¬ 
cients,  can  be  integrated  numerically  through  standard  techniques.  How¬ 
ever,  the  presence  of  the  diffusion  term  makes  integration  from  the  shock 
front  to  the  origin  very  unstable.  The  term  has  a  corresponding  stabiliz¬ 
ing  effect  for  outward  integration.  Unfortunately,  the  natural  specifica¬ 
tion  of  boundary  conditions  at  x=l  calls  for  inward  integration.  In 
principle,  it  is  possible  to  solve  the  equations  by  choosing  parametric 
inner  starting  conditions  and  searching  for  that  set  of  parameters  which 
lead  to  matching  of  the  boundary  conditions  upon  outward  integration.  In 
practice,  the  coupled  nature  of  the  boundary  conditions  in  conjunction 
with  numerical  inaccuracy  makes  this  very  difficult  because  a  nearly 
singular  set  of  conditions  must  be  met. 

Since  the  five  coupled  equations  are  linear  and  homogeneous,  it 
is  possible,  at  least  in  principle,  to  find  a  set  of  five  independent 
solutions  which  can  be  summed  with  with  arbitrary  coefficents  to  form  the 
general  solution.  The  boundary  conditions  are  linear  in  the  dependent 
variables  so,  given  the  general  solution,  the  coefficients  can  be  solved 
for  using  a  simple  system  of  linear  equations.  The  only  problem  is  thus 
in  finding  a  set  of  independent  solutions. 

Independent  solutions  will  give  different  small  x  behavior.  The 
character  of  the  small  x  behavior  of  the  zeroth  order  solutions  is  known; 
this  can  be  used  to  determine  possible  behavior  of  the  perturbed  solu¬ 
tions.  Once  this  is  determined,  numerical  continuations  out  to  the 


boundary  and  appropriate  combinations  to  fit  the  Hugoniot  conditions  along 
with  momentum  conservation  can  be  found.  The  zeroth  order  solutions  have 
the  following  form: 

vi  “  vxx3  +  ... 

Pi  -  p  (1  +  u*2  +  ...) 

(47) 

Pi  “  Px  (1  +  VX4  +  .  .  .  ) 

T i  a  T^(l  -  yX2  +  ...) 

where  the  sub-x  quantities  along  with  y  and  v  are  known  constants.  Using 
the  leading  order  terms  along  with  the  equations  of  motion,  a  set  of  five 

arbitrary  coefficients.  A,  B,  C,  0  and  E,  define  the  small  x  behavior  is 

available  from  independent  solutions  to  the  perturbed  equations: 

gv  =  A  x-1/2  +  B  +  E  X"3 

g^  =  -3/4  A  x-1/2  -  B  +  1/2  E  x~3 

%  1  C  x  +  D  x-2  ♦  E,  (48) 

gp  =  a  Px/5px  B  x  +  apx/2  Ex"2 

=  d  (a  p  /5p  B  -  C/T  )  x  +  p  (-D/T  +  a  p  /2p  E)x2 

pxxx  X  xxxx 

gq  =  c  -  2  d  x-3 

where  t  is  a  nonzero  constant  which  can  be  expressed  in  terms  of  the 

xeroth  order  parameters.  Of  course,  each  of  the  five  independent  solu¬ 
tions,  A-E,  will  produce  some  nonzero  behavior  for  each  of  the  dependent 
variables.  However,  only  the  dominant  or  divergent  terms  are  exhibited 
above.  The  five  independent  solutions  have  been  found  numerically  by 
integration  outward  from  some  small  value  of  x  using  the  small  x  behavior 
from  each  in  turn  as  the  initial  condition.  In  some  cases  it  is  necessary 
to  specify  the  leading  order  term  for  each  of  the  five  dependent  variables 
in  order  to  assure  getting  independent  solutions. 
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There  are  six  constants,  A-E  plus  g^  ,  available  which  are  to  be 
determined  with  four  boundary  conditions  at  x=l,  a  net  zero  momentum  and 
the  requirement  of  a  continuous  temperature  at  the  origin.  The  last 
condition  requires  that  there  be  no  contribution  from  the  D  solution.  It 
is  trivial  to  determine  the  constants  from  the  resulting  linear  conditions 
and  the  solutions  are  shown  in  Figures  5-8  for  the  particular  constants 
indicated  in  the  previous  section.  The  solution  has  no  contribution  from 
the  E  component  and  gives  q^  =  0.69  a2/4H.  This  means  that  the  perturbed 
temperature  vanishes  at  the  origin.  Since  g^.  is  multiplied  by  cos  e  to 
form  the  temperature,  g^  must  vanish  at  x=0,  otherwise  there  would  be  a 
temperature  discontinuity  along  the  axis  through  the  origin.  Such  a  dis¬ 
continuity  would  be  unphysical  since  the  diffusion  term  must  smooth  out 
this  behavior. 


SECTION  4 

VORTICITY  AND  VELOCITY  AT  PRESSURE  EQUILIBRIUM 


The  vorticity  generated  in  the  interior  of  the  early-time  fire¬ 
ball  by  the  passage  of  the  blast  wave  can  result  in  late-time  flow  fields 
which  may  be  important  for  subsequent  dynamics.  Aspects  of  this  behavior 
have  been  given  in  detail  in  our  previous  report11.  We  now  briefly  review 
the  key  features  of  that  report. 


The  vorticity  in  the  flow  field  comes  exclusively  from  the 
exponential  perturbation  of  the  ambient  atmosphere,  that  is  the  first 
order  flow  calculated  in  the  previous  section.  This  vorticity  is 


<o  =  v  X  V  =  (f  /r  +  f  Jr  +  3f  /arlsine  eA 

Q  V  i  (J) 

=  t3/5  u(x)  sine  eA 


(49) 


where  w(x)  contains  the  residual  scaled  dependence  on  position.  At  small 
ranges,  x  «  1,  we  note  that  the  zeroth  order  temperature  is  nearly  con¬ 
stant  according  to  equation  (47).  Since  the  convection  of  the  vorticity 
with  the  flow  is  generally  described  by 


where  S 

vanish. 

m 

in  the  $ 

(a/at  +  v*v)  to/ p  =  <o/o  •  w  +  vT  x  vS 


(50) 
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3/3tff3/5io/ol  =  v/of-o/r  -  Warl'3/5 


(51) 
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The  behavior  of  v  and  p  forces 


o>(x)  =  u)0x-3/2  (52) 

where  <d0  is  a  constant,  if  the  vorticity  convection  is  to  be  satisfied  for 
no  temperature  gradient.  The  small  x  behavior  of  w  is,  in  fact,  exactly 
this  based  on  the  use  of  equation  (47).  Equation  (51)  thus  means  w  can  be 
written  as  a  function  of  r  alone,  independent  of  t.  This  means  that,  at 
small  x,  the  vorticity  is  constant  at  a  fixed  position  given  by  r  and  e 
during  the  strong  blast  wave  phase  of  the  motion. 

Following  the  strong  blast  phase,  a  residual  fireball  remains 
which  is  nearly  isothermal.  The  fireball  expands  until  it  is  nearly  in 
pressure  equilibrium  with  the  surrounding  air.  At  this  point  the  initial 
flow  fields  are  greatly  altered  throuqh  the  doing  of  work  on  expansion; 
however,  there  will  be  residual  flow  associated  with  the  blast  wave  gener¬ 
ated  vorticity  imposed  by  the  initial  gradient  in  the  ambient  atmosphere. 
It  is  this  pressure  equilibrium  flow  field  which  we  wish  to  estimate. 

The  altered  vorticity  distribution  resulting  from  the  nearly 
uniform  expansion  of  the  isothermal  fireball  can  be  found  from  equation 
(50)  by  taking  vT=0  and  using  w  from  equation  (49)  as  the  initial  vortic¬ 
ity.  The  pre-pressure  equilibrium  velocity  field  now  must  be  that  which 
is  dictated  by  fireball  expansion;  we  ignore  the  contribution  of  the  per¬ 
turbed  flow.  Under  this  condition  the  vorticity  equation  can  be  written 
as 


0/Dt(aj/p/r)  =  0 


(53) 


so  that  the  final  m  distribution  can  be  expressed  in  terms  of  the  uniform 
expansion  and  the  intial  distribution.  As  indicated  in  Reference  11  the 
final  vorticity  at  pressure  equilibrium  becomes 


(54) 


wf  =  w0(pf/pi)1/6(a/r)3/2  sine 

where  the  density  ratio,  p^/p.,  is  that  of  the  final  expanded  fireball  to 
that  before  expansion  but  after  the  blast  wave  has  set  up  the  vorticity. 

The  generation  of  velocities  from  the  final  vorticity  field  can 
be  accomplished  by  noting  that  once  pressure  equilibrium  has  been  reached, 
the  subsequent  flow  is  nearly  incompressible.  Consequently,  the  residual 
rotational  velocity  field  can  be  written  in  terms  of  a  vector  potential, 
A,  where 

v  =  v  x  A  (55) 

and  A  is  given  in  terms  of  the  vorticity  as 

A  =  (4tt  )~ 1  /  d3r'  w  ( r ' ) / 1 r-r *  |  (56) 

so  that  given  v  can  be  determined. 

The  estimated  velocity  fields  as  well  as  the  vorticity  are 
smaller  than  suggested  by  the  results  in  Reference  11  by  a  factor  of  about 
ten;  the  previous  calculation  overestimated  <i)0  by  this  amount.  However 
this  just  scales  the  results  so  the  structure  seen  before  was  correct. 
The  calculation  of  Reference  11  differed  from  the  current  one  in  that  it 
assumed  a  zeroth  order  temperature  profile  rather  than  including  the  radi¬ 
ation  diffusion  term.  Examples  of  the  pressure  distribution  velocity 
profiles  found  using  the  radiation  diffusion  term  for  the  40  kT  example  at 
low  altitude  are  shown  in  Figures  9-11  for  three  different  polar  angles. 
The  central  vertical  velocity  is  about  8  meters/second. 
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SECTION  5 
DISCUSSION 


A  similarity  solution  for  the  first  order  perturbation  of  a 
nearly  spherical  strong  blast  wave,  generated  by  a  point  source  in  a 
mildly  stratified  medium,  has  been  found.  The  hydrodynamic  equations,  in 
two  spatial  dimensions,  include  a  radiation  diffusion  approximation  to 
radiation  transport  such  that  it  is  appropriate  for  the  problem  of  low 
altitude  atmospheric  explosions.  The  form  for  radiation  diffusion  is  of 
such  a  character  to  allow  both  the  unperturbed  and  perturbed  solutions  to 
exhibit  a  similarity  form.  The  resulting  cyl indr ical ly  symmetric  numeri¬ 
cal  solution  provides  the  flow  in  both  radial  and  polar  directions  as  well 
as  a  description  of  the  distortion  of  the  nonspherical  shock  front.  It  is 
determined  that  the  stratified  atmosphere  allows  an  increase  in  shock 
speed  in  the  direction  of  decreasing  density  and  to  lowest  order,  the 
additive  correction  factor  can  be  written  as  a  constant  coefficent  times 
a?/4H  »  t4//t>  *  cos  e  where  9  is  the  polar  angle  and  t  is  the  time.  The 
simplified  calculations  by  Kompaneets  gave  a  constant  coefficent  of  one 
while  the  better  calculation  of  Laumbach  and  Probstein,  which  still 
assumed  radial  flow,  gave  the  constant  as  0.75.  The  current  calculation 
gives  0.69  indicating  that  the  inclusion  of  nonradial  flow  tends  to  reduce 
the  distortion  of  the  shock  front. 

The  perturbed  flow  fields  have  nonzero  vorticity,  in  contrast  to 
the  radial  flow  case.  This  vorticity  is  generated  by  the  passage  of  the 
blast  wave  through  the  nonuniform  atmosphere.  Since  the  hot  interior  of 
the  fireball  tends  to  have  a  rather  uniform  temperature,  any  circulation 
generated  before  the  material  passes  inside  will  remain  and  be  convected 
along  with  the  larger  zeroth  order  flow  field.  This  can  provide  a  means 


for  determining  the  flow  fields  in  some  of  the  material  at  times  after  the 
strong  shock  assumption  or  the  perturbation  expansion  cease  to  be  valid. 
The  residual  vorticity  may  be  the  dominant  flow  at  times  the  order  of  one 
second  after  the  fireball  has  ceased  expanding  and  reached  an  approximate 
pressure  equilibrium.  The  vorticity  then  may  play  a  significant  role  in 
the  subsequent  rise  and  torusing  of  the  fireball. 


The  use  of  a  perturbation  expansion  limits  the  value  of  the  cur¬ 
rent  result  to  cases  and  times  for  which  the  shock  radius  is  small  com¬ 
pared  to  a  scale  height.  However,  the  application  of  a  similarity  solu¬ 
tion  qreatly  reduces  the  work  required.  It  might  be  wondered  if  the 
result  here  for  a  perturbation  similarity  solution  can  be  extended  so  as 
to  permit  the  calculation  of  higher  order  terms  to  expand  the  set  of  cases 
which  can  be  studied.  The  answer  is  yes  and  it  can  be  demonstrated  that 
all  orders  in  the  perturbation  expansion  can  be  expressed  in  similarity 
form  so  that  the  solution  is  the  product  of  a  power  of  time  multiplied  by 
a  function  of  the  same  similarity  variable  along  with  a  further  function 
of  polar  angle.  The  resulting  problem  is  tractable  but  it  is  not  as  nice 
as  the  first  order  case  in  that  the  angular  terms  now  contain  contribu¬ 
tions  to  lower  Fourier  orders.  For  example,  the  second  order  terms  supply 
terms  independent  of  angle  as  well  as  terms  proportional  to  cos(2e). 
Furthermore,  the  perturbation  equations  become  even  more  ponderous  with 
increasing  order  as  one  would  expect.  The  appendix  provides  details  of 
the  second  order  perturbation  equations. 
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APPENDIX 

EXTENSION  OF  SIMILARITY  PERTURBATION  TECHNIQUE 


A.l  GENERAL  FORM. 

As  illustrated  in  the  main  body  of  the  note,  it  is  possible  to 
solve  the  first  order  pertubation  problem  of  a  strong  blast  wave  in  a 
mildly  stratified  atmosphere  by  finding  the  zeroth  order  solution,  which 
is  known  to  permit  a  similarity  solution,  and  using  it  as  the  basis  for  a 
perturbation  to  the  lowest  order  in  the  stratification  parameter.  The 
remarkable  feature  of  these  first  order  equations  is  that  they  also  admit 
a  similarity  solution  in  the  same  similarity  variable  as  used  for  zeroth 
order.  In  this  appendix  we  shall  demonstrate  that  this  result  continues 
to  hold  for  all  orders  of  the  perturbation  expansion  and  the  character  of 
the  required  functions  are  provided  at  arbitrary  order;  the  resulting 
equations  for  the  second  order  perturbation  will  explicitly  written  down. 
For  this  purpose  we  shall  ignore  the  complications  associated  with  the 
radiation  diffusion  term.  However,  the  same  line  of  reasoning  holds  in 
that  case  as  well  if  the  form  of  the  opacity  is  chosen  so  as  to  permit  the 
similarity  solution  at  zeroth  order. 

The  hydrodynamic  equations  outside  the  shock  are 


where  viscosity  is  ignored  and  an  ideal  gas  equation  of  state 


p  =  PR*T  (58) 

is  assumed. 

For  a  strong  shock  the  boundary  conditions  which  must  be  met  at 
the  (unknown)  location  of  the  shock  front  are 

p(s)  _  y+1 
Pa  y-1 

vi(s)  =  f.  2n(s)..)1/2  (59) 

(y+1)p3 

v,(s)  0 

3G  _1 _  _  /(y+l)p(s)  >1/2 

9t  | vG|  2Pa 

where  G(r,t)  =  0  defines  the  shock  location,  s,  and  pa  is  the  ambient 

a 

density. 

The  parallel  and  perpendicular  directions  are  relative  to  the 
local  shock  front. 

For  a  uniform  atmosphere  with  zero  external  ambient  pressure 
(which  defines  the  strong  shock  condition)  the  solution  must  be  express- 
able  in  similarity  form  since  there  are  no  length  scales  associated  with 
the  problem.  It  can  easily  be  shown  that  the  solution  can  be  written  as 

p  =  p0(x)f6/5 

D  =  Po (x)  (60) 

v  =  v0(x)t-3/5 
r 
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and 


G  =  at2/5  -  r 


where  a  is  dictated  by  the  total  energy  in  the  problem. 

The  similarity  variable  x  is  the  ratio  of  the  radius,  r,  to  the 
shock  radius,  at2/5,  and  the  time  dependence  is  separable.  By  symmetry, 
there  is  no  angular  dependence. 

Consider  the  problem  imposed  by  a  pertubation  in  the  ambient 
density  so  that  it  is  no  longer  uniform  but  is  a  function  of  only  one 
cartesian  dimension.  For  a  realistic  atmosphere,  the  dependence  is 
usual ly  taken  as 

pa  =  Pb  exp( -z/H)  (61) 

where  the  scale  height  H  completely  characterizes  the  variation.  (For  the 
following  discussion,  any  mild  analytic  variation  of  density  which  can  be 
expressed  in  terms  of  a  single  length  parameter  is  equivalent:  the  result 
is  not  unique  to  exponential  density  variation.)  The  scale  height  enters 
the  problem  exclusively  through  the  appearence  of  ambient  density  in 
boundary  conditions. 

The  expansion  of  the  ambient  density  is 


p  =  P  (1  -  r...cose  +  1  rji_cos9)  -  ...) 

a  b  H  2  H 


(62) 


where  the  natural  expansion  parameter  which  appears  is  r/H  and  so  the 
usual  small  perturbation  restriction  suggests  that  the  truncated  results 
can  be  expected  to  valid  only  if  the  radius  of  the  shock  is  much  less  than 
the  scale  height. 


The  angular  dependence  of  the  perturbation  terms  in  the  density 
are  simply  factorable  powers  of  cos(e).  Consequently,  since  the  perturba¬ 
tion  expansion  hydrodynamic  equations,  and  boundary  conditions,  involve  at 
most  multiplication  and  differentiation,  it  will  always  be  possible  to 
express  the  solutions  any  order  in  terms  of  a  truncated  Fourier  series; 
the  number  of  Fourier  terms  must  increase  with  increasing  perturbation 
order . 

Consider  a  candidate  series  for  the  perturbation  series  which 
has  the  similarity  form  with  individual  terms  showing  separable  angular 
and  time  dependence: 


p  =  p0(x)  +  Pl(x,e)tPl  +  p2(x,e)tp2 


p  =  p0(x)t-6/5  +  p^x.eHP1  +  p2(x,e)tP2  + 


v  r  =  v0(x)t-3/5  +  v^x.ejt71  +  v2(x,9)tv2  + 
Oi(x,e)t91  +  e2(x,e)te2  + 


(63) 


ve  = 


and  an  auxiliary  function,  G,  giving  the  position  of  the  shock  front 


G  =  at2/5  +  G ! ( 0 ) tG1  +  G2 (0 ) t'3^  +  ...  -r 


G2 


(64) 


where  r  is  the  shock  position  at  polar  angle  e  for  G(r,e,t)=0.  (Note  that 
vG  defines  the  normal  to  the  shock  front.)  The  similarity  variable  is 
x  =  r/R(t)  where  R(t)=at2/5  is  the  zeroth  order  shock  position.  We  now 
ask  if  such  a  similarity  form  can  provide  a  solution  to  the  equations 
together  with  boundary  conditions. 
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Consider  that  the  proposed  solution  is  inserted  into  the  full 
equations  and  the  resulting  equations  are  expanded  and  written  as  a  sum  of 
terms,  each  beinq  a  particular  order  in  the  perturbation  series,  each  as 
coefficent  of  powers  of  the  expansion  parameter,  1/H.  If  our  proposed 
similarity  condition  is  to  be  met,  the  individual  orders  must  have  the 
same  factorable  time  dependence  in  each  of  its  contributions  from  all 
appropriate  orders  of  the  dependent  variables.  However,  the  different 
orders  may  not,  in  fact,  cannot,  have  the  same  time  dependence. 

In  order  for  this  to  be  true  certain  restrictions  are  required 
based  on  the  character  of  the  original  equations.  For  example,  the  con¬ 
tinuity  equation  will  require  that  the  nth  order  density  contribution 
less  one  power  of  time  must  have  the  time  power  as  the  product  of  the 
ith  or(jer  in  density  and  the  {n-i)th  order  of  velocity  aside  from  the 
effects  of  the  divergence  operation.  It  is  easily  seen  that  the  diver¬ 
gence  operation  can  be  expressed  as  an  operation  in  the  x  and  e  variables 
along  with  a  division  by  R(t)  which  is  independent  of  order.  Consequent¬ 
ly,  since  the  i=n  case  is  proper  (-1,  for  s/at  =  -3/5,  for  v0  plus  -2/5, 
for  1/R( t)  --  all  of  which  is  independent  of  the  nth  order  density 
power),  it  remains  only  to  assure  that  this  holds  for  n*i.  It  is  evident 
that  this  requires  that  the  increments  in  the  powers  of  time  for  succes¬ 
sive  orders  of  both  density  and  velocity  must  be  equal  and  equal  for  all 
orders.  Inspection  of  the  other  hydrodynamic  equations  indicates  this 
must  be  true  for  all  the  dependent  variables.  Note  that  this  behavior 
forces  the  product  of  an  arbitrary  number  of  such  series  to  hold  the  same 
property  -  that  of  uniform  time  power  increments  between  orders. 

Note  that  we  have  not  yet  imposed  any  conditions  concerning  the 
boundary  behavior,  which  must  include  ambient  density  variation.  We  have 
only  determined  an  allowable  form  for  a  perturbation  expansion  which  meets 
the  similarity  condition.  It  remains  to  be  seen  that  this  is  consistent 
with  the  perturbed  density  and  boundary  conditions.  The  expansion  of  the 
ambient  density  can  be  written  as 


which  suggests  that  the  desired  power  of  time  increment  must  be  2/5.  How¬ 
ever,  it  must  be  remembered  that  the  shock  conditions  are  imposed  at  the 
actual  shock  front  whose  perturbed  location  (not  x=l)  is  to  be  determined 
as  a  part  of  the  solution.  Generally  the  evaluation  of  an  nth  order 
hydrodynamic  variable,  say  the  pressure,  at  the  shock  will  have  the  form: 

p  +  y  f.tfi  9p-P-i  +  y  y  f,  t^f.t^j  _ g-k-j  +  y  y  y  ... 

n  i  1  ar  kj  k  J  ar2  (66) 

all  evaluated  at  x=l.  Since  the  a/ar  supplies  a  extra  t-2/5  dependence, 
all  terms  at  nth  order  will  have  the  same  power  of  time  if  the  increment 
in  the  powers  of  the  expansion  of  G  is  2/5  just  as  for  the  ambient  dens¬ 
ity. 

Consequently,  the  variables  evaluated  at  the  shock  front  all 
have  the  desired  behavior  of  associating  a  single  power  of  time  with  a 
particular  order  while  incrementing  that  power  uniformly  by  2/5  with  each 
order.  Since  the  leading  order  has  the  appropriate  power  for  all  of  the 
boundary  conditions,  all  orders  must  also.  Therefore  the  postulated  form 
for  perturbation  similarity  solutions  exists  for  the  problem  at  hand.  The 
dependence  of  the  dependent  variables  at  arbitrary  nth  order  can  now  be 
written  as 


P  ' 


Pn(x,0)t2n/5 

2n/5-6/5 


Pn  ( x , a )  t 


vn(x,a)t 


2n/5-3/5 


Gn(e)t 


2n/5+2/5 


(67) 


$ 


In  fact,  as  indicated  before,  the  angular  dependence  for  each  term  is 
separable  in  the  above  expressions  and  the  nt(l  order  dependence  can  be 
expressed  as  a  truncated  Fourier  series  terminating  at  cos(ne)  or  sin(ne) 
dependina  on  the  parity  of  the  variable. 

A. 2  SPECIFIC  EQUATIONS  AT  2nd  ORDER. 


VI 


The  form  and  solution  for  the  first  order  perturbation  terms 
have  been  given  in  the  main  text  for  a  specific  case.  At  higher  orders 
the  expansions  becomes  increasingly  complex  and  it  is  doubtful  if  there  is 
utility  in  forminq  the  general  term.  However,  it  is  illustrative  to  go 
one  step  further  and  explicitly  provide  the  second  order  equations  since 
this  brings  out  the  essential  features  of  the  angular  dependence.  To  be 
consistent  with  the  main  text  the  expansion  is  taken  to  be 


p(r,t)  =  p0(x)  +  gp(x)  f2/5  cose  +  hp 
p(r,t)  =  p0(x)’6/5  +  gp(x)  t'u/5  cose  +  hp 

vr(r,t)  =  v0(x)-3/5  +  gr(x)  t-i/5  cos0  + 


vQ(r,t)  - 


9e(x)  t'1/5  cose  +  h0 


g(r,t)  =  R ( t)  +  gR(x)  t"4/5cose  +  hR(t,e)  -  r 


■'.>1 


*<*:•)* 


‘if*:*. 


The  second  order  perturbed  equations  are  thus 


i.  (r2hf)+  i  (SM|  h,  *  m)}  *  cosVp[L  ^  (r2fr)  *  l  f.l 

h  [1_  i_  (r2v)  *  v  it  *  cos28f.  It  +  h -  slnia.  f  f  =  o 
p  r2  3p  ar  r  ar  rar  r  6  p 

Yp[l_  ( r2hr  1+  I  he  +  ^M)l  +  YC0S2efp[-L.  i-  (r2fr)  +  1  fQ] 


r2  ar  r  r  sine 


yML  J-  C2v)  *  V  ♦  cos2efr  L»  *  ^  f,fp  ■  0 


P  r2  3p 


v  — L  +  cos2e  f  — £  +  h  —  -  faiffl+fr)  S1'n20 
ar  r  ar  r  3r  e  e  r 

.  3d  hn  f2  -  ,  fp  3fp  i  3hp 

.  1  ap  f-JL  +  cos 2 9 ]  +  cos2e - -  -  i  — 

d  ar  P  p2  p2  3r  p  ar 


aft  vhQ  sine  cose 


f  ah„  ,  afQ  vh 

_1  sine  cose  +  v  — i  +  sine  cose  f — 2.  +  — 2.  + 
r  ar  rar  r 


1  ahD  fn  fD 

i _ —  -  J-  S  sine  cose 

rp  ae  rp  p 


The  second  order  boundary  conditions  are 


8  f 

3^  +  cose  ^  f2  +  cos20  f  p  +  h  =  xti  p.  [-  !r  cos2e  +  iL  cos2e) 
R  ar  2  ar2  R  R  ar  P  y-1  bl  H  ?H2 


-sine  coseff  +  11  f  )  +  v  ~ B.  =  -R f h  +  fDcose  sine)  -  fDfnsine  cose 
r  3r  k  39  9  ar  K  K  9 


+  1  ,L_v  f2  cos2e  +  — B.  cos2e  +  h  +  sj_n2Q.  f  f, 


R  o  2  R 


r  r  e  R 


=  (  I  tI^h.+  J^f2COS29  +  I-^f0COS29 

VI  0b  v  P  ar  R  2P  ar2  R  Pa»  ' 

+  £flS-2.l  r R  in  +i  )f  +  cos2e  +  sin2!  f2  +  cos.,9  Li  r  + 

H  par  R  2H2  R2  R  H  p  p 


2  a  r,,  3P  ,  .  f 


((i  f  )2  ♦  [^)2  h-  (i)2)} 

n  ^  m  •'  n  H 


8  par 


— -B-  =  x±L  fiv  hD  +  1  liv.  f2  cos2e  +  fD  -i-L  cos2e  +  h  +  sin2e  f  f  ) 
at  2  ■  ar  R  2  ar2  R  R  ar  r  R  9  R 


(70) 


The  similarity  form  for  the  solution  provides  a  h  time  dependence  which 
allows  the  equations  to  be  simplified.  Furthermore,  the  angular  depend¬ 
ence  can  be  shown  to  allow  only  the  following  form 

h  =  t1/5  H  sine  cose 

n  e 

h  =  t4/5  (H  cos2e  +  J  ) 

0  D  P 

I  hp  =  t*2/5  (Hp  cos2e  +  Jp)  (71) 

hr  =  t1/,5(Hr  sine  cose  +  Jr) 
hR  =  t6/ 5 ( hr  sine  cose  +  JR)  . 


Using  both  the  time  and  angular  separable  dependence,  the  equations  can 
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